In [ ]:
import stim
import pymatching
import numpy as np
import sinter
from typing import List
import matplotlib.pyplot as plt
In [ ]:
def generate_circuit(rounds=1, before_round_data_depolarization= 0.1, before_measure_flip_probability=0.2, after_reset_flip_probability=0.3, after_clifford_depolarization=0.4):
return f"""# surface_code circuit exercise.
# rounds: {rounds}
# before_round_data_depolarization: {before_round_data_depolarization}
# before_measure_flip_probability: {before_measure_flip_probability}
# after_reset_flip_probability: {after_reset_flip_probability}
# after_clifford_depolarization: {after_clifford_depolarization}
# layout:
# L00 X01 L02 X03 L04 X05 L06
# Z10 d11 Z12 d13 Z14 d15 Z16
# d20 X21 d22 X22 d24 X25 d26
# Z30 d31 Z32 d33 Z34 d35 Z36
# d40 X41 d42 X42 d44 X45 d46
# Z50 d51 Z52 d53 Z54 d55 Z56
# d60 X61 d62 X62 d64 X65 d66
# Legend:
# d# = data qubit
# L# = data qubit with logical observable crossing
# X# = measurement qubit (X stabilizer)
# Z# = measurement qubit (Z stabilizer)
QUBIT_COORDS(0, 0) 1
QUBIT_COORDS(0, 1) 2
QUBIT_COORDS(0, 2) 3
QUBIT_COORDS(0, 3) 4
QUBIT_COORDS(0, 4) 5
QUBIT_COORDS(0, 5) 6
QUBIT_COORDS(0, 6) 7
QUBIT_COORDS(1, 0) 8
QUBIT_COORDS(1, 1) 9
QUBIT_COORDS(1, 2) 10
QUBIT_COORDS(1, 3) 11
QUBIT_COORDS(1, 4) 12
QUBIT_COORDS(1, 5) 13
QUBIT_COORDS(1, 6) 14
QUBIT_COORDS(2, 0) 15
QUBIT_COORDS(2, 1) 16
QUBIT_COORDS(2, 2) 17
QUBIT_COORDS(2, 3) 18
QUBIT_COORDS(2, 4) 19
QUBIT_COORDS(2, 5) 20
QUBIT_COORDS(2, 6) 21
QUBIT_COORDS(3, 0) 22
QUBIT_COORDS(3, 1) 23
QUBIT_COORDS(3, 2) 24
QUBIT_COORDS(3, 3) 25
QUBIT_COORDS(3, 4) 26
QUBIT_COORDS(3, 5) 27
QUBIT_COORDS(3, 6) 28
QUBIT_COORDS(4, 0) 29
QUBIT_COORDS(4, 1) 30
QUBIT_COORDS(4, 2) 31
QUBIT_COORDS(4, 3) 32
QUBIT_COORDS(4, 4) 33
QUBIT_COORDS(4, 5) 34
QUBIT_COORDS(4, 6) 35
QUBIT_COORDS(5, 0) 36
QUBIT_COORDS(5, 1) 37
QUBIT_COORDS(5, 2) 38
QUBIT_COORDS(5, 3) 39
QUBIT_COORDS(5, 4) 40
QUBIT_COORDS(5, 5) 41
QUBIT_COORDS(5, 6) 42
QUBIT_COORDS(6, 0) 43
QUBIT_COORDS(6, 1) 44
QUBIT_COORDS(6, 2) 45
QUBIT_COORDS(6, 3) 46
QUBIT_COORDS(6, 4) 47
QUBIT_COORDS(6, 5) 48
QUBIT_COORDS(6, 6) 49
R 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
X_ERROR({after_reset_flip_probability}) 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
R 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
X_ERROR({after_reset_flip_probability}) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
TICK
DEPOLARIZE1({before_round_data_depolarization}) 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
TICK
# Plaquette Syndromes
H 2 4 6 16 18 20 30 32 34 44 46 48
DEPOLARIZE1({after_clifford_depolarization}) 2 4 6 16 18 20 30 32 34 44 46 48
TICK
# Plaquette Syndromes triplets
CX 1 2 3 4 5 6 43 44 45 46 47 48
DEPOLARIZE2({after_clifford_depolarization}) 1 2 3 4 5 6 43 44 45 46 47 48
TICK
CX 9 2 11 4 13 6 37 44 39 46 41 48
DEPOLARIZE2({after_clifford_depolarization}) 9 2 11 4 13 6 37 44 39 46 41 48
TICK
CX 3 2 5 4 7 6 45 44 47 46 49 48
DEPOLARIZE2({after_clifford_depolarization}) 3 2 5 4 7 6 45 44 47 46 49 48
TICK
# Plaquette Syndromes quads
CX 9 16 11 18 13 20 23 30 25 32 27 34
DEPOLARIZE2({after_clifford_depolarization}) 9 16 11 18 13 20 23 30 25 32 27 34
TICK
CX 15 16 17 18 19 20 29 30 31 32 33 34
DEPOLARIZE2({after_clifford_depolarization}) 15 16 17 18 19 20 29 30 31 32 33 34
TICK
CX 17 16 19 18 21 20 31 30 33 32 35 34
DEPOLARIZE2({after_clifford_depolarization}) 17 16 19 18 21 20 31 30 33 32 35 34
TICK
CX 23 16 25 18 27 20 37 30 39 32 41 34
DEPOLARIZE2({after_clifford_depolarization}) 23 16 25 18 34 20 37 30 39 32 41 34
TICK
# Vertex Syndromes triplets
CX 1 8 15 22 29 36 7 14 21 28 35 42
DEPOLARIZE2({after_clifford_depolarization}) 1 8 15 22 29 36 7 14 21 28 35 42
TICK
CX 15 8 29 22 43 36 21 14 35 28 49 42
DEPOLARIZE2({after_clifford_depolarization}) 1 8 15 22 29 36 7 14 21 28 35 42
TICK
CX 9 8 23 22 37 36 13 14 27 28 41 42
DEPOLARIZE2({after_clifford_depolarization}) 1 8 15 22 29 36 7 14 21 28 35 42
TICK
# Vertex Syndromes quads
CX 9 10 11 12 23 24 25 26 37 38 39 40
DEPOLARIZE2({after_clifford_depolarization}) 9 10 11 12 23 24 25 26 37 38 39 40
TICK
CX 3 10 5 12 17 24 19 26 31 38 33 40
DEPOLARIZE2({after_clifford_depolarization}) 3 10 5 12 17 24 19 26 31 38 33 40
TICK
CX 11 10 13 12 25 24 27 26 39 38 41 40
DEPOLARIZE2({after_clifford_depolarization}) 11 10 13 12 25 24 27 26 39 38 41 40
TICK
CX 17 10 19 12 31 24 33 26 45 38 47 40
DEPOLARIZE2({after_clifford_depolarization}) 17 10 19 12 31 24 33 26 45 38 47 40
TICK
# Plaquette Syndromes
H 2 4 6 16 18 20 30 32 34 44 46 48
DEPOLARIZE1({after_clifford_depolarization}) 2 4 6 16 18 20 30 32 34 44 46 48
TICK
# Measure measurement qubits
X_ERROR({before_measure_flip_probability}) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
MR 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
X_ERROR({after_reset_flip_probability}) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
# Detect errors using the vertex Syndromes
DETECTOR(1, 0, 0) rec[-21]
DETECTOR(3, 0, 0) rec[-14]
DETECTOR(5, 0, 0) rec[-7]
DETECTOR(1, 2, 0) rec[-20]
DETECTOR(3, 2, 0) rec[-13]
DETECTOR(5, 2, 0) rec[-6]
DETECTOR(1, 4, 0) rec[-19]
DETECTOR(3, 4, 0) rec[-12]
DETECTOR(3, 4, 0) rec[-5]
DETECTOR(1, 6, 0) rec[-18]
DETECTOR(3, 6, 0) rec[-11]
DETECTOR(5, 6, 0) rec[-4]
REPEAT {rounds}""" + "{" + f"""
TICK
DEPOLARIZE1({before_round_data_depolarization}) 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
TICK
# Plaquette Syndromes
H 2 4 6 16 18 20 30 32 34 44 46 48
DEPOLARIZE1({after_clifford_depolarization}) 2 4 6 16 18 20 30 32 34 44 46 48
TICK
# Plaquette Syndromes triplets
CX 1 2 3 4 5 6 43 44 45 46 47 48
DEPOLARIZE2({after_clifford_depolarization}) 1 2 3 4 5 6 43 44 45 46 47 48
TICK
CX 9 2 11 4 13 6 37 44 39 46 41 48
DEPOLARIZE2({after_clifford_depolarization}) 9 2 11 4 13 6 37 44 39 46 41 48
TICK
CX 3 2 5 4 7 6 45 44 47 46 49 48
DEPOLARIZE2({after_clifford_depolarization}) 3 2 5 4 7 6 45 44 47 46 49 48
TICK
# Plaquette Syndromes quads
CX 9 16 11 18 13 20 23 30 25 32 27 34
DEPOLARIZE2({after_clifford_depolarization}) 9 16 11 18 13 20 23 30 25 32 27 34
TICK
CX 15 16 17 18 19 20 29 30 31 32 33 34
DEPOLARIZE2({after_clifford_depolarization}) 15 16 17 18 19 20 29 30 31 32 33 34
TICK
CX 17 16 19 18 21 20 31 30 33 32 35 34
DEPOLARIZE2({after_clifford_depolarization}) 17 16 19 18 21 20 31 30 33 32 35 34
TICK
CX 23 16 25 18 27 20 37 30 39 32 41 34
DEPOLARIZE2({after_clifford_depolarization}) 23 16 25 18 34 20 37 30 39 32 41 34
TICK
# Vertex Syndromes triplets
CX 1 8 15 22 29 36 7 14 21 28 35 42
DEPOLARIZE2({after_clifford_depolarization}) 1 8 15 22 29 36 7 14 21 28 35 42
TICK
CX 15 8 29 22 43 36 21 14 35 28 49 42
DEPOLARIZE2({after_clifford_depolarization}) 1 8 15 22 29 36 7 14 21 28 35 42
TICK
CX 9 8 23 22 37 36 13 14 27 28 41 42
DEPOLARIZE2({after_clifford_depolarization}) 1 8 15 22 29 36 7 14 21 28 35 42
TICK
# Vertex Syndromes quads
CX 9 10 11 12 23 24 25 26 37 38 39 40
DEPOLARIZE2({after_clifford_depolarization}) 9 10 11 12 23 24 25 26 37 38 39 40
TICK
CX 3 10 5 12 17 24 19 26 31 38 33 40
DEPOLARIZE2({after_clifford_depolarization}) 3 10 5 12 17 24 19 26 31 38 33 40
TICK
CX 11 10 13 12 25 24 27 26 39 38 41 40
DEPOLARIZE2({after_clifford_depolarization}) 11 10 13 12 25 24 27 26 39 38 41 40
TICK
CX 17 10 19 12 31 24 33 26 45 38 47 40
DEPOLARIZE2({after_clifford_depolarization}) 17 10 19 12 31 24 33 26 45 38 47 40
TICK
# Plaquette Syndromes
H 2 4 6 16 18 20 30 32 34 44 46 48
DEPOLARIZE1({after_clifford_depolarization}) 2 4 6 16 18 20 30 32 34 44 46 48
TICK
# Measure measurement qubits
X_ERROR({before_measure_flip_probability}) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
MR 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
X_ERROR({after_reset_flip_probability}) 2 4 6 8 10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 44 46 48
# Used for adding time?
SHIFT_COORDS(0, 0, 1)
# Detecting errors by comparing to previous shot
DETECTOR(6, 5, 0) rec[-1] rec[-25]
DETECTOR(6, 3, 0) rec[-2] rec[-26]
DETECTOR(6, 1, 0) rec[-3] rec[-27]
DETECTOR(5, 6, 0) rec[-4] rec[-28]
DETECTOR(5, 4, 0) rec[-5] rec[-29]
DETECTOR(5, 2, 0) rec[-6] rec[-30]
DETECTOR(5, 0, 0) rec[-7] rec[-31]
DETECTOR(4, 5, 0) rec[-8] rec[-32]
DETECTOR(4, 3, 0) rec[-9] rec[-33]
DETECTOR(4, 1, 0) rec[-10] rec[-34]
DETECTOR(3, 6, 0) rec[-11] rec[-35]
DETECTOR(3, 4, 0) rec[-12] rec[-36]
DETECTOR(3, 2, 0) rec[-13] rec[-37]
DETECTOR(3, 0, 0) rec[-14] rec[-38]
DETECTOR(2, 5, 0) rec[-15] rec[-39]
DETECTOR(2, 3, 0) rec[-16] rec[-40]
DETECTOR(2, 1, 0) rec[-17] rec[-41]
DETECTOR(1, 6, 0) rec[-18] rec[-42]
DETECTOR(1, 4, 0) rec[-19] rec[-43]
DETECTOR(1, 2, 0) rec[-20] rec[-44]
DETECTOR(1, 0, 0) rec[-21] rec[-45]
DETECTOR(0, 5, 0) rec[-22] rec[-46]
DETECTOR(0, 3, 0) rec[-23] rec[-47]
DETECTOR(0, 1, 0) rec[-24] rec[-48]
"""+"}"+f"""
# Final readout of the 25 data qubits in Z basis as we are using vertex syndromes
X_ERROR({before_measure_flip_probability}) 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
M 1 3 5 7 9 11 13 15 17 19 21 23 25 27 29 31 33 35 37 39 41 43 45 47 49
DETECTOR(1, 0, 1) rec[-18] rec[-21] rec[-25] rec[-46]
DETECTOR(3, 0, 1) rec[-11] rec[-14] rec[-18] rec[-39]
DETECTOR(5, 0, 1) rec[-4] rec[-7] rec[-11] rec[-32]
DETECTOR(1, 2, 1) rec[-17] rec[-20] rec[-21] rec[-24] rec[-45]
DETECTOR(3, 2, 1) rec[-10] rec[-13] rec[-14] rec[-17] rec[-38]
DETECTOR(5, 2, 1) rec[-3] rec[-6] rec[-7] rec[-10] rec[-31]
DETECTOR(1, 4, 1) rec[-16] rec[-19] rec[-20] rec[-23] rec[-44]
DETECTOR(3, 4, 1) rec[-9] rec[-12] rec[-13] rec[-16] rec[-37]
DETECTOR(5, 4, 1) rec[-2] rec[-5] rec[-6] rec[-9] rec[-30]
DETECTOR(1, 6, 1) rec[-15] rec[-19] rec[-22] rec[-43]
DETECTOR(3, 6, 1) rec[-8] rec[-12] rec[-15] rec[-36]
DETECTOR(5, 6, 1) rec[-1] rec[-5] rec[-8] rec[-29]
OBSERVABLE_INCLUDE(0) rec[-22] rec[-23] rec[-24] rec[-25]
"""
In [ ]:
circuit = stim.Circuit(generate_circuit())
In [ ]:
circuit.without_noise().diagram('timeslice-svg')
In [ ]:
def count_logical_errors(circuit: stim.Circuit, num_shots: int) -> int:
# Sample the circuit.
sampler = circuit.compile_detector_sampler()
sample =sampler.sample(num_shots, separate_observables=True)
detection_events, observable_flips = sample
detection_events = np.array(detection_events, order='C')
# Configure a decoder using the circuit.
detector_error_model = circuit.detector_error_model(decompose_errors=True)
matcher = pymatching.Matching.from_detector_error_model(detector_error_model)
# Run the decoder.
predictions = []
for i in range(num_shots):
predictions.append(matcher.decode(detection_events[i]))
predictions=np.array(predictions)
# Count the mistakes.
num_errors = 0
for shot in range(num_shots):
actual_for_shot = observable_flips[shot]
predicted_for_shot = predictions[shot]
if not np.array_equal(actual_for_shot, predicted_for_shot):
num_errors += 1
return num_errors
In [ ]:
circuit = stim.Circuit(generate_circuit(rounds=4))
num_shots = 100_000
num_logical_errors = count_logical_errors(circuit, num_shots)
print("there were", num_logical_errors, "wrong predictions (logical errors) out of", num_shots, "shots")
In [ ]:
dem = circuit.detector_error_model()
print(repr(dem))
In [ ]:
sampler = circuit.compile_sampler()
sample = sampler.sample(shots=1)
one_sample = sample[0]
for k in range(0, len(one_sample), 8):
timeslice = one_sample[k:k+8]
print("".join("1" if e else "_" for e in timeslice))
In [ ]:
circuit.without_noise().diagram(
"detslice-with-ops-svg",
tick=range(0, 50),
)
In [ ]:
circuit.diagram("matchgraph-3d")
In [ ]:
surface_code_tasks = [
sinter.Task(
circuit = stim.Circuit(
generate_circuit(
rounds=3*4,
after_clifford_depolarization=noise,
after_reset_flip_probability=noise,
before_measure_flip_probability=noise,
before_round_data_depolarization=noise,
)
),
json_metadata={'d': d, 'r': d * 3, 'p': noise},
)
for d in [2]
for noise in [0.0001,0.005,0.001,0.01,0.02,0.03,0.06,0.09,0.1,0.2,0.3]
]
collected_surface_code_stats: List[sinter.TaskStats] = sinter.collect(
num_workers=4,
tasks=surface_code_tasks,
decoders=['pymatching'],
max_shots=1_000_000,
max_errors=5_000,
print_progress=True,
)
In [ ]:
fig, ax = plt.subplots(1, 1)
sinter.plot_error_rate(
ax=ax,
stats=collected_surface_code_stats,
x_func=lambda stat: stat.json_metadata['p'],
group_func=lambda stat: stat.json_metadata['d'],
failure_units_per_shot_func=lambda stat: stat.json_metadata['r'],
)
ax.loglog()
ax.set_title("Surface Code Error Rates per Round under Circuit Noise")
ax.set_xlabel("Phyical Error Rate")
ax.set_ylabel("Logical Error Rate per Round")
ax.grid(which='major')
ax.legend()
fig.set_dpi(120) # Show it bigger
In [ ]:
from qccd import *
# Create architecture
arch = QCCDArch()
ionspacing = 5
trapspacing = 30
cooling_colour = 'red'
qubit_colour = 'lightblue'
traps = []
for i in range(3):
ions = [QubitIon(qubit_colour, 'Q'), QubitIon(qubit_colour, 'Q'), QubitIon(qubit_colour, 'Q')]
trap = arch.addManipulationTrap(x=i *trapspacing, y=0, ions=ions, color='grey', spacing=ionspacing, isHorizontal=True)
traps.append(trap)
for t1, t2 in zip(traps[:-1], traps[1:]):
arch.addEdge(t1, t2)
arch.refreshGraph()
fig, ax = plt.subplots()
# Display architecture
arch.display(fig,ax, showLabels=False)
In [ ]:
from qccd import *
# Create architecture
arch = QCCDArch()
ionspacing = 1
trapspacing = 6
qubit_colour = 'red'
qubit_colour = 'lightblue'
junction_colour = 'orange'
traps = []
ManipulationTrap.DEFAULT_SPACING = trapspacing
arch.SIZING = 0.5
ions11 = [QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q'), QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q')]
trap11 = arch.addManipulationTrap(x=0, y=0, ions=ions11, color='grey', spacing=ionspacing, isHorizontal=False)
ions12 = [QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q'), QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q')]
trap12 = arch.addManipulationTrap(x=trapspacing, y=0, ions=ions12, color='grey', spacing=ionspacing, isHorizontal=False)
junctionL = arch.addJunction(x=0, y=trapspacing, color=junction_colour)
junctionR = arch.addJunction(x=trapspacing, y=trapspacing, color=junction_colour)
ions21 = [QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q'), QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q')]
trap21 = arch.addManipulationTrap(x=0, y=2*trapspacing, ions=ions21, color='grey', spacing=ionspacing, isHorizontal=False)
ions22 = [QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q'), QubitIon(qubit_colour, 'Q'), CoolingIon(qubit_colour, 'Q')]
trap22 = arch.addManipulationTrap(x=trapspacing, y=2*trapspacing, ions=ions22, color='grey', spacing=ionspacing, isHorizontal=False)
crossing11 = arch.addEdge(trap11, junctionL)
arch.addEdge(trap12, junctionR)
crossing21 = arch.addEdge(trap21, junctionL)
arch.addEdge(trap22, junctionR)
arch.addEdge(junctionL, junctionR)
arch.refreshGraph()
# Display architecture
ops = (
Split.physicalOperation(trap11, crossing11),
Move.physicalOperation(crossing11),
JunctionCrossing.physicalOperation(junctionL, crossing11),
JunctionCrossing.physicalOperation(junctionL, crossing21),
Move.physicalOperation(crossing21),
Merge.physicalOperation(trap21, crossing21),
GateSwap.physicalOperation(trap=trap21, ion1=ions21[0], ion2=ions21[2]),
TwoQubitMSGate.physicalOperation(ion1=ions22[0], ion2=ions22[2],trap=trap22),
OneQubitGate.physicalOperation(ion=ions12[0], trap=trap12),
Measurement.physicalOperation(ion=ions12[0], trap=trap12),
QubitReset.physicalOperation(ion=ions12[0],trap=trap12))
parallelOps = paralleliseOperations(ops)
ncols = 2
fig, ax = plt.subplots()
arch.display(fig, ax, showLabels=False)
fig.set_size_inches(arch.WINDOW_SIZE[0]*0.5, arch.WINDOW_SIZE[1]*0.5)
fig,axs = plt.subplots(int(len(parallelOps)/2), 2)
axs = axs.flatten()
for i, (ax, parallelOp)in enumerate(zip(axs, parallelOps)):
parallelOp.run()
arch.refreshGraph()
print(parallelOp.label+f' fidelity {parallelOp.fidelity()} operation time {parallelOp.operationTime()}')
arch.display(fig=fig, ax=ax, title=parallelOp.label, operation=parallelOp, showLabels=False)
fig.set_size_inches(arch.WINDOW_SIZE[0]*0.5, arch.WINDOW_SIZE[1]*int(len(parallelOps)/2)*0.5)
In [ ]:
from qccd import *
noise = 1e-3
d=4
circuit = QCCDCircuit.generated(
"surface_code:rotated_memory_z",
rounds=1,
distance=d,
after_clifford_depolarization=noise,
after_reset_flip_probability=noise,
before_measure_flip_probability=noise,
before_round_data_depolarization=noise,
)
rows = 4
cols = 4
capacity = 4
arch, instructions = circuit.processCircuit(rows=rows, cols=cols, trapCapacity=capacity, isHorizontal=(rows==1), clustering=1)
for i in instructions:
print(i.label+str([ion.idx for ion in i.ions]))
capacities = [1,2]
rows = 6
cols = 6
label='part'
fig, axs = plt.subplots(1, len(capacities))
for ax, capacity in zip(axs, capacities):
arch1, instructions = circuit.processCircuit(rows=rows, cols=cols, trapCapacity=capacity, isHorizontal=(rows==1), clustering=2)
arch1.ION_SIZE = 2000
for idx, (ion, pos) in circuit._ionMapping.items():
ion.set(ion.idx, pos[0], pos[1])
arch1.refreshGraph()
edgesDups = []
edgesPos = []
ionsInvolved = set()
score = 0
for op in instructions:
if not isinstance(op, TwoQubitMSGate):
continue
if not ionsInvolved.isdisjoint(op.ions):
score+=1
ionsInvolved=set()
edgesDups.append((op.ions, score))
edgesPos.append((op.ionsActedIdxs, score))
ionsInvolved=ionsInvolved.union(op.ions)
scores = [score for ((ion1,ion2), score) in edgesDups]
color_from_score = {s: (5+i*i)*2 for i, s in enumerate( sorted(list(set(scores)), reverse=True))}
arch1._manipulationTraps.append(([(ion1.idx, ion2.idx) for ((ion1,ion2), score) in edgesDups], [color_from_score[score] for ((ion1,ion2), score) in edgesDups]))
arch1.display(fig, ax, showLabels=False, showEdges=False, show_junction=False)
ax.set_title(f"{label} {capacities}")
circuit.without_noise().diagram('timeline-svg')
QUBIT_RESET[2] QUBIT_RESET[3] QUBIT_RESET[15] QUBIT_RESET[16] QUBIT_RESET[26] QUBIT_RESET[27] QUBIT_RESET[39] QUBIT_RESET[6] QUBIT_RESET[19] QUBIT_RESET[30] QUBIT_RESET[42] QUBIT_RESET[9] QUBIT_RESET[33] QUBIT_RESET[34] QUBIT_RESET[45] QUBIT_RESET[12] QUBIT_RESET[1] QUBIT_RESET[14] QUBIT_RESET[25] QUBIT_RESET[38] QUBIT_RESET[5] QUBIT_RESET[18] QUBIT_RESET[29] QUBIT_RESET[41] QUBIT_RESET[8] QUBIT_RESET[21] QUBIT_RESET[32] QUBIT_RESET[44] QUBIT_RESET[11] QUBIT_RESET[23] QUBIT_RESET[36] RY[1] RX[1] RY[14] RX[14] RY[38] RX[38] RY[29] RX[29] RY[8] RX[8] RY[44] RX[44] RY[23] RX[23] RY[36] RX[36] RY[1] RX[1] RX[3] TWO_QUBIT_MS_GATE[1, 3] RY[1] RY[29] RX[29] RX[30] TWO_QUBIT_MS_GATE[29, 30] RY[29] RY[38] RX[38] RX[39] TWO_QUBIT_MS_GATE[38, 39] RY[38] RY[44] RX[44] RX[45] TWO_QUBIT_MS_GATE[44, 45] RY[44] RY[14] RX[14] RX[16] TWO_QUBIT_MS_GATE[14, 16] RY[14] RY[8] RX[8] RX[9] TWO_QUBIT_MS_GATE[8, 9] RY[8] RY[19] RX[19] RX[18] TWO_QUBIT_MS_GATE[19, 18] RY[19] RY[27] RX[27] RX[25] TWO_QUBIT_MS_GATE[27, 25] RY[27] RY[34] RX[34] RX[32] TWO_QUBIT_MS_GATE[34, 32] RY[34] RY[42] RX[42] RX[41] TWO_QUBIT_MS_GATE[42, 41] RY[42] RY[6] RX[6] RX[5] TWO_QUBIT_MS_GATE[6, 5] RY[6] RY[12] RX[12] RX[11] TWO_QUBIT_MS_GATE[12, 11] RY[12] RY[1] RX[1] RX[2] TWO_QUBIT_MS_GATE[1, 2] RY[1] RY[29] RX[29] RX[19] TWO_QUBIT_MS_GATE[29, 19] RY[29] RY[38] RX[38] RX[27] TWO_QUBIT_MS_GATE[38, 27] RY[38] RY[44] RX[44] RX[34] TWO_QUBIT_MS_GATE[44, 34] RY[44] RY[14] RX[14] RX[15] TWO_QUBIT_MS_GATE[14, 15] RY[14] RY[8] RX[8] RX[42] TWO_QUBIT_MS_GATE[8, 42] RY[8] RY[26] RX[26] RX[18] TWO_QUBIT_MS_GATE[26, 18] RY[26] RY[3] RX[3] RX[25] TWO_QUBIT_MS_GATE[3, 25] RY[3] RY[30] RX[30] RX[32] TWO_QUBIT_MS_GATE[30, 32] RY[30] RY[39] RX[39] RX[41] TWO_QUBIT_MS_GATE[39, 41] RY[39] RY[16] RX[16] RX[5] TWO_QUBIT_MS_GATE[16, 5] RY[16] RY[9] RX[9] RX[11] TWO_QUBIT_MS_GATE[9, 11] RY[9] RY[29] RX[29] RX[27] TWO_QUBIT_MS_GATE[29, 27] RY[29] RY[23] RX[23] RX[34] TWO_QUBIT_MS_GATE[23, 34] RY[23] RY[38] RX[38] RX[15] TWO_QUBIT_MS_GATE[38, 15] RY[38] RY[44] RX[44] RX[42] TWO_QUBIT_MS_GATE[44, 42] RY[44] RY[8] RX[8] RX[6] TWO_QUBIT_MS_GATE[8, 6] RY[8] RY[36] RX[36] RX[12] TWO_QUBIT_MS_GATE[36, 12] RY[36] RY[26] RX[26] RX[25] TWO_QUBIT_MS_GATE[26, 25] RY[26] RY[33] RX[33] RX[32] TWO_QUBIT_MS_GATE[33, 32] RY[33] RY[30] RX[30] RX[41] TWO_QUBIT_MS_GATE[30, 41] RY[30] RY[39] RX[39] RX[5] TWO_QUBIT_MS_GATE[39, 5] RY[39] RY[45] RX[45] RX[11] TWO_QUBIT_MS_GATE[45, 11] RY[45] RY[9] RX[9] RX[21] TWO_QUBIT_MS_GATE[9, 21] RY[9] RY[29] RX[29] RX[26] TWO_QUBIT_MS_GATE[29, 26] RY[29] RY[23] RX[23] RX[33] TWO_QUBIT_MS_GATE[23, 33] RY[23] RY[38] RX[38] RX[3] TWO_QUBIT_MS_GATE[38, 3] RY[38] RY[44] RX[44] RX[30] TWO_QUBIT_MS_GATE[44, 30] RY[44] RY[8] RX[8] RX[39] TWO_QUBIT_MS_GATE[8, 39] RY[8] RY[36] RX[36] RX[45] TWO_QUBIT_MS_GATE[36, 45] RY[36] RY[2] RX[2] RX[25] TWO_QUBIT_MS_GATE[2, 25] RY[2] RY[19] RX[19] RX[32] TWO_QUBIT_MS_GATE[19, 32] RY[19] RY[27] RX[27] RX[41] TWO_QUBIT_MS_GATE[27, 41] RY[27] RY[15] RX[15] RX[5] TWO_QUBIT_MS_GATE[15, 5] RY[15] RY[42] RX[42] RX[11] TWO_QUBIT_MS_GATE[42, 11] RY[42] RY[6] RX[6] RX[21] TWO_QUBIT_MS_GATE[6, 21] RY[6] RY[1] RX[1] RY[14] RX[14] RY[38] RX[38] RY[29] RX[29] RY[8] RX[8] RY[44] RX[44] RY[23] RX[23] RY[36] RX[36] MEASUREMENT[1] MEASUREMENT[14] MEASUREMENT[25] MEASUREMENT[38] MEASUREMENT[5] MEASUREMENT[18] MEASUREMENT[29] MEASUREMENT[41] MEASUREMENT[8] MEASUREMENT[21] MEASUREMENT[32] MEASUREMENT[44] MEASUREMENT[11] MEASUREMENT[23] MEASUREMENT[36] QUBIT_RESET[1] QUBIT_RESET[14] QUBIT_RESET[25] QUBIT_RESET[38] QUBIT_RESET[5] QUBIT_RESET[18] QUBIT_RESET[29] QUBIT_RESET[41] QUBIT_RESET[8] QUBIT_RESET[21] QUBIT_RESET[32] QUBIT_RESET[44] QUBIT_RESET[11] QUBIT_RESET[23] QUBIT_RESET[36] MEASUREMENT[2] MEASUREMENT[3] MEASUREMENT[15] MEASUREMENT[16] MEASUREMENT[26] MEASUREMENT[27] MEASUREMENT[39] MEASUREMENT[6] MEASUREMENT[19] MEASUREMENT[30] MEASUREMENT[42] MEASUREMENT[9] MEASUREMENT[33] MEASUREMENT[34] MEASUREMENT[45] MEASUREMENT[12]
Out[ ]:
In [ ]:
from qccd import *
noise = 1e-3
d=4
circuit = QCCDCircuit.generated(
"surface_code:rotated_memory_z",
rounds=1,
distance=d,
after_clifford_depolarization=noise,
after_reset_flip_probability=noise,
before_measure_flip_probability=noise,
before_round_data_depolarization=noise,
)
rows = 4
cols = 4
capacity = 4
arch, instructions = circuit.processCircuit(rows=rows, cols=cols, trapCapacity=capacity, isHorizontal=(rows==1), clustering=1)
for i in instructions:
print(i.label+str([ion.idx for ion in i.ions]))
capacities = [2,3,4,5,6,7]
rows = 6
cols = 6
for i, label in enumerate(('basic', 'pos', 'part')):
fig, axs = plt.subplots(1, len(capacities))
for ax, capacity in zip(axs, capacities):
arch1, instructions = circuit.processCircuit(rows=rows, cols=cols, trapCapacity=capacity, isHorizontal=(rows==1), clustering=i)
arch1.ION_SIZE = 2000
for idx, (ion, pos) in circuit._ionMapping.items():
ion.set(ion.idx, pos[0], pos[1])
arch1.refreshGraph()
edgesDups = []
edgesPos = []
ionsInvolved = set()
score = 0
for op in instructions:
if not isinstance(op, TwoQubitMSGate):
continue
if not ionsInvolved.isdisjoint(op.ions):
score+=1
ionsInvolved=set()
edgesDups.append((op.ions, score))
edgesPos.append((op.ionsActedIdxs, score))
ionsInvolved=ionsInvolved.union(op.ions)
scores = [score for ((ion1,ion2), score) in edgesDups]
color_from_score = {s: (5+i*i)*0 for i, s in enumerate( sorted(list(set(scores)), reverse=True))}
arch1._manipulationTraps.append(([(ion1.idx, ion2.idx) for ((ion1,ion2), score) in edgesDups], [color_from_score[score] for ((ion1,ion2), score) in edgesDups]))
arch1.display(fig, ax, showLabels=False, showEdges=False, show_junction=False)
ax.set_title(f"{label} {capacities}")
circuit.without_noise().diagram('timeline-svg')
QUBIT_RESET[2] QUBIT_RESET[3] QUBIT_RESET[15] QUBIT_RESET[16] QUBIT_RESET[26] QUBIT_RESET[27] QUBIT_RESET[39] QUBIT_RESET[6] QUBIT_RESET[19] QUBIT_RESET[30] QUBIT_RESET[42] QUBIT_RESET[9] QUBIT_RESET[33] QUBIT_RESET[34] QUBIT_RESET[45] QUBIT_RESET[12] QUBIT_RESET[1] QUBIT_RESET[14] QUBIT_RESET[25] QUBIT_RESET[38] QUBIT_RESET[5] QUBIT_RESET[18] QUBIT_RESET[29] QUBIT_RESET[41] QUBIT_RESET[8] QUBIT_RESET[21] QUBIT_RESET[32] QUBIT_RESET[44] QUBIT_RESET[11] QUBIT_RESET[23] QUBIT_RESET[36] RY[1] RX[1] RY[14] RX[14] RY[38] RX[38] RY[29] RX[29] RY[8] RX[8] RY[44] RX[44] RY[23] RX[23] RY[36] RX[36] RY[1] RX[1] RX[3] TWO_QUBIT_MS_GATE[1, 3] RY[1] RY[29] RX[29] RX[30] TWO_QUBIT_MS_GATE[29, 30] RY[29] RY[38] RX[38] RX[39] TWO_QUBIT_MS_GATE[38, 39] RY[38] RY[44] RX[44] RX[45] TWO_QUBIT_MS_GATE[44, 45] RY[44] RY[14] RX[14] RX[16] TWO_QUBIT_MS_GATE[14, 16] RY[14] RY[8] RX[8] RX[9] TWO_QUBIT_MS_GATE[8, 9] RY[8] RY[19] RX[19] RX[18] TWO_QUBIT_MS_GATE[19, 18] RY[19] RY[27] RX[27] RX[25] TWO_QUBIT_MS_GATE[27, 25] RY[27] RY[34] RX[34] RX[32] TWO_QUBIT_MS_GATE[34, 32] RY[34] RY[42] RX[42] RX[41] TWO_QUBIT_MS_GATE[42, 41] RY[42] RY[6] RX[6] RX[5] TWO_QUBIT_MS_GATE[6, 5] RY[6] RY[12] RX[12] RX[11] TWO_QUBIT_MS_GATE[12, 11] RY[12] RY[1] RX[1] RX[2] TWO_QUBIT_MS_GATE[1, 2] RY[1] RY[29] RX[29] RX[19] TWO_QUBIT_MS_GATE[29, 19] RY[29] RY[38] RX[38] RX[27] TWO_QUBIT_MS_GATE[38, 27] RY[38] RY[44] RX[44] RX[34] TWO_QUBIT_MS_GATE[44, 34] RY[44] RY[14] RX[14] RX[15] TWO_QUBIT_MS_GATE[14, 15] RY[14] RY[8] RX[8] RX[42] TWO_QUBIT_MS_GATE[8, 42] RY[8] RY[26] RX[26] RX[18] TWO_QUBIT_MS_GATE[26, 18] RY[26] RY[3] RX[3] RX[25] TWO_QUBIT_MS_GATE[3, 25] RY[3] RY[30] RX[30] RX[32] TWO_QUBIT_MS_GATE[30, 32] RY[30] RY[39] RX[39] RX[41] TWO_QUBIT_MS_GATE[39, 41] RY[39] RY[16] RX[16] RX[5] TWO_QUBIT_MS_GATE[16, 5] RY[16] RY[9] RX[9] RX[11] TWO_QUBIT_MS_GATE[9, 11] RY[9] RY[29] RX[29] RX[27] TWO_QUBIT_MS_GATE[29, 27] RY[29] RY[23] RX[23] RX[34] TWO_QUBIT_MS_GATE[23, 34] RY[23] RY[38] RX[38] RX[15] TWO_QUBIT_MS_GATE[38, 15] RY[38] RY[44] RX[44] RX[42] TWO_QUBIT_MS_GATE[44, 42] RY[44] RY[8] RX[8] RX[6] TWO_QUBIT_MS_GATE[8, 6] RY[8] RY[36] RX[36] RX[12] TWO_QUBIT_MS_GATE[36, 12] RY[36] RY[26] RX[26] RX[25] TWO_QUBIT_MS_GATE[26, 25] RY[26] RY[33] RX[33] RX[32] TWO_QUBIT_MS_GATE[33, 32] RY[33] RY[30] RX[30] RX[41] TWO_QUBIT_MS_GATE[30, 41] RY[30] RY[39] RX[39] RX[5] TWO_QUBIT_MS_GATE[39, 5] RY[39] RY[45] RX[45] RX[11] TWO_QUBIT_MS_GATE[45, 11] RY[45] RY[9] RX[9] RX[21] TWO_QUBIT_MS_GATE[9, 21] RY[9] RY[29] RX[29] RX[26] TWO_QUBIT_MS_GATE[29, 26] RY[29] RY[23] RX[23] RX[33] TWO_QUBIT_MS_GATE[23, 33] RY[23] RY[38] RX[38] RX[3] TWO_QUBIT_MS_GATE[38, 3] RY[38] RY[44] RX[44] RX[30] TWO_QUBIT_MS_GATE[44, 30] RY[44] RY[8] RX[8] RX[39] TWO_QUBIT_MS_GATE[8, 39] RY[8] RY[36] RX[36] RX[45] TWO_QUBIT_MS_GATE[36, 45] RY[36] RY[2] RX[2] RX[25] TWO_QUBIT_MS_GATE[2, 25] RY[2] RY[19] RX[19] RX[32] TWO_QUBIT_MS_GATE[19, 32] RY[19] RY[27] RX[27] RX[41] TWO_QUBIT_MS_GATE[27, 41] RY[27] RY[15] RX[15] RX[5] TWO_QUBIT_MS_GATE[15, 5] RY[15] RY[42] RX[42] RX[11] TWO_QUBIT_MS_GATE[42, 11] RY[42] RY[6] RX[6] RX[21] TWO_QUBIT_MS_GATE[6, 21] RY[6] RY[1] RX[1] RY[14] RX[14] RY[38] RX[38] RY[29] RX[29] RY[8] RX[8] RY[44] RX[44] RY[23] RX[23] RY[36] RX[36] MEASUREMENT[1] MEASUREMENT[14] MEASUREMENT[25] MEASUREMENT[38] MEASUREMENT[5] MEASUREMENT[18] MEASUREMENT[29] MEASUREMENT[41] MEASUREMENT[8] MEASUREMENT[21] MEASUREMENT[32] MEASUREMENT[44] MEASUREMENT[11] MEASUREMENT[23] MEASUREMENT[36] QUBIT_RESET[1] QUBIT_RESET[14] QUBIT_RESET[25] QUBIT_RESET[38] QUBIT_RESET[5] QUBIT_RESET[18] QUBIT_RESET[29] QUBIT_RESET[41] QUBIT_RESET[8] QUBIT_RESET[21] QUBIT_RESET[32] QUBIT_RESET[44] QUBIT_RESET[11] QUBIT_RESET[23] QUBIT_RESET[36] MEASUREMENT[2] MEASUREMENT[3] MEASUREMENT[15] MEASUREMENT[16] MEASUREMENT[26] MEASUREMENT[27] MEASUREMENT[39] MEASUREMENT[6] MEASUREMENT[19] MEASUREMENT[30] MEASUREMENT[42] MEASUREMENT[9] MEASUREMENT[33] MEASUREMENT[34] MEASUREMENT[45] MEASUREMENT[12]
Out[ ]:
In [ ]:
from qccd_circuit import *
noise=0.01
circuit = QCCDCircuit.generated(
"surface_code:rotated_memory_z",
rounds=1,
distance=4,
after_clifford_depolarization=noise,
after_reset_flip_probability=noise,
before_measure_flip_probability=noise,
before_round_data_depolarization=noise,
)
circuit.without_noise().diagram('timeslice-svg')
Out[ ]:
In [ ]:
arch, instructions = circuit.processCircuitAugmentedGrid(rows=3, cols=3, trapCapacity=7)
oldPositions = {}
for idx, (ion, pos) in circuit._ionMapping.items():
oldPositions[idx] = ion.pos
ion.set(ion.idx, pos[0], pos[1], parent=ion.parent)
fig,ax =plt.subplots()
arch1 = arch
arch1.refreshGraph()
edgesDups = []
edgesPos = []
ionsInvolved = set()
score = 0
for op in instructions:
if not isinstance(op, TwoQubitMSGate):
continue
if not ionsInvolved.isdisjoint(op.ions):
score+=1
ionsInvolved=set()
edgesDups.append((op.ions, score))
edgesPos.append((op.ionsActedIdxs, score))
ionsInvolved=ionsInvolved.union(op.ions)
scores = [score for ((ion1,ion2), score) in edgesDups]
color_from_score = {s: (5+i*i)*0 for i, s in enumerate( sorted(list(set(scores)), reverse=True))}
arch1._manipulationTraps.append(([(ion1.idx, ion2.idx) for ((ion1,ion2), score) in edgesDups], [color_from_score[score] for ((ion1,ion2), score) in edgesDups]))
arch1.display(fig, ax, showLabels=False, showEdges=False, show_junction=False)
fig.set_size_inches(arch.WINDOW_SIZE[0]*0.6, arch.WINDOW_SIZE[1]*0.6)
arch1._manipulationTraps = arch1._manipulationTraps[:-1]
for idx, (ion, pos) in circuit._ionMapping.items():
ion.set(ion.idx, oldPositions[idx][0], oldPositions[idx][1], parent=ion.parent)
arch.refreshGraph()
fig,ax =plt.subplots()
arch.display(fig, ax, showLabels=False)
fig.set_size_inches(arch.WINDOW_SIZE[0]*0.8, arch.WINDOW_SIZE[1]*0.8)
allOps, _ = arch.processOperationsViaForwarding(instructions)
parallelOpsMap = paralleliseOperations(allOps)
parallelOps = list(dict(parallelOpsMap).values())
# arch = circuit.resetArch()
# arch.refreshGraph()
# axsPerFig = 2
# figs = []
# axs = []
# numFigs = int(len(parallelOps)/(axsPerFig))+int((len(parallelOps)%(axsPerFig))>0)
# for _ in range(numFigs):
# fig, axs_ = plt.subplots(1, axsPerFig, sharex=False, sharey=False)
# figs.extend([fig]*axsPerFig)
# axs.extend(axs_)
# figs = figs[:len(parallelOps)]
# axs = axs[:len(parallelOps)]
# for fig, ax, parallelOp in zip(figs, axs, parallelOps):
# arch.display(fig, ax, operation=parallelOp, runOps=True, showLabels=False)
print(f"total number of qubit operations: {len(instructions)}")
print(f"total number of operations: {len(allOps)}")
print(f"time for operations: {max(parallelOpsMap.keys())}")
total number of qubit operations: 349 total number of operations: 610 time for operations: 0.02321800000000003
In [ ]:
from qccd_circuit import *
noise=0.01
circuit = QCCDCircuit.generated(
"surface_code:rotated_memory_z",
rounds=1,
distance=4,
after_clifford_depolarization=noise,
after_reset_flip_probability=noise,
before_measure_flip_probability=noise,
before_round_data_depolarization=noise,
)
arch, instructions = circuit.processCircuitAugmentedGrid(rows=5, cols=5, trapCapacity=1)
oldPositions = {}
for idx, (ion, pos) in circuit._ionMapping.items():
oldPositions[idx] = ion.pos
ion.set(ion.idx, pos[0], pos[1], parent=ion.parent)
fig,ax =plt.subplots()
arch1 = arch
arch1.refreshGraph()
edgesDups = []
edgesPos = []
ionsInvolved = set()
score = 0
for op in instructions:
if not isinstance(op, TwoQubitMSGate):
continue
if not ionsInvolved.isdisjoint(op.ions):
score+=1
ionsInvolved=set()
edgesDups.append((op.ions, score))
edgesPos.append((op.ionsActedIdxs, score))
ionsInvolved=ionsInvolved.union(op.ions)
scores = [score for ((ion1,ion2), score) in edgesDups]
color_from_score = {s: (5+i*i)*0 for i, s in enumerate( sorted(list(set(scores)), reverse=True))}
arch1._manipulationTraps.append(([(ion1.idx, ion2.idx) for ((ion1,ion2), score) in edgesDups], [color_from_score[score] for ((ion1,ion2), score) in edgesDups]))
arch1.display(fig, ax, showLabels=False, showEdges=False, show_junction=False)
fig.set_size_inches(arch.WINDOW_SIZE[0]*0.5, arch.WINDOW_SIZE[1]*0.5)
arch1._manipulationTraps = arch1._manipulationTraps[:-1]
for idx, (ion, pos) in circuit._ionMapping.items():
ion.set(ion.idx, oldPositions[idx][0], oldPositions[idx][1], parent=ion.parent)
arch.refreshGraph()
fig,ax =plt.subplots()
arch.display(fig, ax, title='map complete', showLabels=False)
fig.set_size_inches(arch.WINDOW_SIZE[0]*0.5, arch.WINDOW_SIZE[1]*0.5)
allOps, barriers = arch.processOperationsViaForwarding(instructions)
parallelOpsMap = paralleliseOperationsWithBarriers(allOps, barriers)
parallelOps = list(dict(parallelOpsMap).values())
circuit.simulate(allOps)
arch = circuit.resetArch()
arch.refreshGraph()
print(f"total number of qubit operations: {len(instructions)}")
print(f"total number of operations: {len(allOps)}")
print(f"time for operations: {max(parallelOpsMap.keys())}")
total number of qubit operations: 349 total number of operations: 955 time for operations: 0.012392
In [ ]:
print(len(paralleliseOperationsSimple(allOps)))
In [ ]:
parallelOps = paralleliseOperationsSimple(allOps)
arch = circuit.resetArch()
arch.refreshGraph()
axsPerFig = 2
figs = []
axs = []
numFigs = int(len(parallelOps)/(axsPerFig))+int((len(parallelOps)%(axsPerFig))>0)
for _ in range(numFigs):
fig, axs_ = plt.subplots(1, axsPerFig, sharex=False, sharey=False)
figs.extend([fig]*axsPerFig)
axs.extend(axs_)
figs = figs[:len(parallelOps)]
axs = axs[:len(parallelOps)]
for fig, ax, parallelOp in zip(figs, axs, parallelOps):
arch.display(fig, ax, operation=parallelOp, runOps=True, showLabels=False)
In [ ]:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
%matplotlib inline
from IPython.display import HTML, clear_output, display
arch = circuit.resetArch()
arch.refreshGraph()
fig,ax=plt.subplots()
parallelOpsSimple = paralleliseOperationsSimple(allOps)
def update(frame, _arch):
ax.clear()
# Clear the output for dynamic display in Jupyter notebook
clear_output(wait=True)
if frame>0:
op: ParallelOperation = parallelOpsSimple[frame-1]
title = f"Operation: {op.label}"
_arch.display(fig, ax, title, operation=op, runOps=True, showLabels=False)
else:
_arch.display(fig, ax, showLabels=False)
return display(fig)
import time
time.sleep(10)
ani = FuncAnimation(fig, lambda frame: update(frame, arch), frames=len(parallelOpsSimple)+1, repeat=False)
ani.save('animation.mp4', writer='ffmpeg', fps=10, dpi=100)
# HTML(ani.to_jshtml())
In [ ]:
import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation
%matplotlib inline
from IPython.display import HTML, clear_output, display
arch = circuit.resetArch()
arch.refreshGraph()
fig,ax=plt.subplots()
parallelOpsTimes = sorted(parallelOpsMap.keys())
def update(frame, _arch):
ax.clear()
# Clear the output for dynamic display in Jupyter notebook
clear_output(wait=True)
if frame>0:
op: ParallelOperation = parallelOps[frame-1]
title = f"Operation {round(parallelOpsTimes[frame-1],6)}: {op.label}"
_arch.display(fig, ax, title, operation=op, runOps=True, showLabels=False)
else:
_arch.display(fig, ax, showLabels=False)
return display(fig)
import time
time.sleep(10)
ani = FuncAnimation(fig, lambda frame: update(frame, arch), frames=len(parallelOps)+1, repeat=False)
ani.save('animation.mp4', writer='ffmpeg', fps=10, dpi=100)
HTML(ani.to_jshtml())
In [ ]:
# gate improvement that is expected in the next few years is about 10X or maximum 100X, not more than that
gate_improvements = [0.025, 0.05, 0.1, 0.25, 0.5, 1.0, 5.0, 10.0, 12.5, 15.0, 20.0, 25.0, 35.0, 50.0, 100.0]
distances = [3,5,7,9,11]
# capacities refer to maximum number of qubits we place in a single trap for the initial QCCD configuration
# true maximum trap capacity, i.e. the maximum number of ions each trap can hold in any instance, is 5 times this initial configuration maximum value
capacities = [1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,30,60,120]
data = {k: {"Forwarding": {d: {c: None for c in capacities} for d in distances}, "Routing": {d: {c: None for c in capacities} for d in distances}} for k in ("PhysicalErrorRates", "Operations", "ElapsedTime", "MeanConcurrency", "QubitOperations", "LogicalErrorRates")}
In [ ]:
for k1 in data.keys():
for k2 in data[k1].keys():
for k3 in data[k1][k2].keys():
print(k1, k2, k3, data[k1][k2][k3])
In [ ]:
data = ...
In [ ]:
with open('saved_Data.txt', 'w') as f:
f.write(str(data))
In [ ]:
from qccd_circuit import *
num_shots = 1_000_000
data: Dict[str, Dict[str, Dict[int, Dict[int, Any]]]]
import concurrent.futures
import logging
import os
num_cores = os.cpu_count()
num_processes = num_cores
logger = logging.getLogger()
logger.setLevel(logging.INFO)
handler = logging.FileHandler("process_log.txt")
formatter = logging.Formatter('%(processName)s - %(levelname)s - %(message)s')
handler.setFormatter(formatter)
logger.addHandler(handler)
with concurrent.futures.ProcessPoolExecutor(max_workers=num_processes) as executor:
futures = []
for d in distances:
for c in capacities:
logger.info(f"Submitting task for distance {d} and capacity {c}.")
futures.append(executor.submit(process_circuit, d, c, gate_improvements, num_shots))
for future in concurrent.futures.as_completed(futures):
if future.exception() is not None:
logger.error(f"Exception Occured for future", exc_info=future.exception())
continue
try:
result = future.result()
d = result["Distance"]
c = result["Capacity"]
for label in result["ElapsedTime"]:
data["ElapsedTime"][label][d][c] = result["ElapsedTime"][label]
data["Operations"][label][d][c] = result["Operations"][label]
data["MeanConcurrency"][label][d][c] = result["MeanConcurrency"][label]
data["QubitOperations"][label][d][c] = result["QubitOperations"][label]
data["LogicalErrorRates"][label][d][c] = result["LogicalErrorRates"][label]
data["PhysicalErrorRates"][label][d][c] = result["PhysicalErrorRates"][label]
logger.info(f"Processing result for distance {result['Distance']}, capacity {result['Capacity']}.")
except Exception as e:
logger.error(f"Exception Occured for future", exc_info=True)
continue
In [ ]:
for k1 in data.keys():
for k2 in data[k1].keys():
for k3 in data[k1][k2].keys():
data[k1][k2][k3] = [data[k1][k2][k3][c] for c in capacities]
data: Dict[str, Dict[str, Dict[int, Dict[int, Any]]]]
In [ ]:
import matplotlib as mpl
import matplotlib.colors as mplcolors
In [ ]:
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, 3))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, ax1 = plt.subplots()
iLabel, label = 0, "Forwarding"
ax1.plot(distances, [data["Operations"][label][d][0] for d in distances], c=colors[0], linestyle=("dashed" if iLabel==1 else "solid"), label=f"Total Operations for QEC shot")
ax1.plot(distances, [int(data["Operations"][label][d][0]/data["MeanConcurrency"][label][d][0]) for d in distances], c=colors[1], linestyle=("dashed" if iLabel==1 else "solid"), label=f"Mean Circuit Depth for QEC shot")
ax1.scatter(distances, [data["Operations"][label][d][0] for d in distances], c=colors[0], marker=("v" if iLabel==1 else "o"))
ax1.scatter(distances, [int(data["Operations"][label][d][0]/data["MeanConcurrency"][label][d][0]) for d in distances], c=colors[1], marker=("v" if iLabel==1 else "o"))
ax1.legend()
ax1.set_title('Importance of Effective Scheduling')
ax1.set_xticks(distances)
ax1.set_xlabel('Distances')
ax1.set_ylabel('Num. Ops.')
ax1.grid(which='major')
In [ ]:
_distances = [3,5,7,9,11]
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, 3))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, ax1 = plt.subplots()
iLabel, label = 0, "Forwarding"
ax1.plot(_distances, [data["Operations"][label][d][0] for d in _distances], c=colors[0], linestyle=("dashed" if iLabel==1 else "solid"), label=f"Total Number of Operations")
ax1.plot(_distances, [data["QubitOperations"][label][d][0] for d in _distances], c=colors[1], linestyle=("dashed" if iLabel==1 else "solid"), label=f"Number of Gates")
ax1.scatter(_distances, [data["Operations"][label][d][0] for d in _distances], c=colors[0], marker=("v" if iLabel==1 else "o"))
ax1.scatter(_distances, [data["QubitOperations"][label][d][0] for d in _distances], c=colors[1], marker=("v" if iLabel==1 else "o"))
ax1.legend()
ax1.set_title('Importance of Effective Mapping')
ax1.set_xlabel('Distances')
ax1.set_xticks(_distances)
ax1.set_ylabel('Operations')
ax1.grid(which='major')
In [ ]:
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(distances)+1))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, ax1 = plt.subplots()
iLabel, label = 0, "Forwarding"
for j, d in enumerate(distances):
ax1.plot(capacities, [data["Operations"][label][d][i] for i in range(len(capacities))], c=colors[j], linestyle=("dashed" if iLabel==1 else "solid"), label=f"d={d}")
ax1.scatter(capacities, [data["Operations"][label][d][i] for i in range(len(capacities))], c=colors[j], marker=("v" if iLabel==1 else "o"))
ax1.legend()
ax1.set_title('Total Number of Operations in One Shot of a QEC Cycle')
ax1.set_xlabel('Trap Capacity')
ax1.set_ylabel('Num. Ops.')
ax1.grid(which='major')
fig1, ax1 = plt.subplots()
iLabel, label = 0, "Forwarding"
for j, d in enumerate([distances[0], distances[-1]]):
ax1.plot(capacities, [data["MeanConcurrency"][label][d][i] for i in range(len(capacities))], c=colors[j], linestyle=("dashed" if iLabel==1 else "solid"), label=f"d={d}")
ax1.scatter(capacities, [data["MeanConcurrency"][label][d][i] for i in range(len(capacities))], c=colors[j], marker=("v" if iLabel==1 else "o"))
ax1.legend()
ax1.set_title('Concurrent Operations in a Timeslice of a QEC Cycle')
ax1.set_xlabel('Trap Capacity')
ax1.set_ylabel('Num. Ops.')
ax1.grid(which='major')
In [ ]:
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(capacities)+1))
colors = [mplcolors.to_hex(c) for c in colors]
fig5, ax5 = plt.subplots()
for i, c in enumerate(capacities):
if i%4==0:
ax5.plot(distances, [round(data["MeanConcurrency"][label][d][i]/(d**2-1), 3) for d in distances], c=colors[i], linestyle=("dashed" if iLabel==1 else "solid"), label=f"c={c}")
ax5.scatter(distances, [round(data["MeanConcurrency"][label][d][i]/(d**2-1), 3) for d in distances], c=colors[i], marker=("v" if iLabel==1 else "o"))
ax5.legend()
ax5.set_title('Concurrency vs Capacity')
ax5.set_xticks(distances)
ax5.set_xlabel('Code Distance')
ax5.set_ylabel('Concurrent Operations per Plaquette')
ax5.grid(which='major')
fig4, ax4 = plt.subplots()
for i, c in enumerate(capacities):
if i%4==0:
ax4.plot(distances, [data["ElapsedTime"][label][d][i] for d in distances], c=colors[i], linestyle=("dashed" if iLabel==1 else "solid"),label=f"c={c}")
ax4.scatter(distances, [data["ElapsedTime"][label][d][i] for d in distances], c=colors[i], marker=("v" if iLabel==1 else "o"))
ax4.legend()
ax4.set_title('Mean Time for One Shot of a QEC Cycle')
ax4.set_xlabel('Distance')
ax4.set_ylabel('Time Elapsed')
ax4.grid(which='major')
In [ ]:
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(distances)+1))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, axs1 = plt.subplots(int(len(capacities)/2),2)
axs1=axs1.flatten()
fig1.set_size_inches( 2*fig4.get_size_inches()[0],fig4.get_size_inches()[1]*len(capacities)/1.4)
thresholds = []
thresholdIdxs = []
for j, (c, ax1) in enumerate(zip(capacities, axs1)):
threshold = gate_improvements[0]
thresholdIdx = 0
for idx, improv in enumerate(gate_improvements[1:]):
isThreshold = True
for d1, d2 in zip(distances[:-1], distances[1:]):
if data["LogicalErrorRates"]["Forwarding"][d1][j][idx+1]<data["LogicalErrorRates"]["Forwarding"][d2][j][idx+1]:
isThreshold = False
if isThreshold:
threshold = (gate_improvements[idx-1]+gate_improvements[idx])/2
thresholdIdx = idx
break
thresholds.append(threshold)
thresholdIdxs.append(thresholdIdx)
for i,d in enumerate(distances):
ax1.plot(gate_improvements, data["LogicalErrorRates"]["Forwarding"][d][j], linestyle='dashed', c=colors[i])
ax1.scatter(gate_improvements, data["LogicalErrorRates"]["Forwarding"][d][j], c=colors[i], label=f"d={d}")
ax1.vlines([threshold], ymin=[1e-7], ymax=[1e-0], colors=[colors[-1]], label="Threshold Point")
ax1.legend()
ax1.set_title(f'Threshold={threshold} for capacity={c}')
ax1.set_xlabel('Gate Improvement')
ax1.set_ylabel('Logical Error Rate')
ax1.set_xlim(1e-2, 1000)
ax1.set_ylim(1e-7, 1e-0)
ax1.set_xticks(gate_improvements)
ax1.loglog()
ax1.grid(which='major')
In [ ]:
print(thresholds)
print(thresholdIdxs)
In [ ]:
thresholds[capacities.index(120)] = 7.5
thresholds[capacities.index(60)] = 7.5
thresholds[capacities.index(30)] = 7.5
thresholds[capacities.index(15)] = 7.5
thresholds[capacities.index(14)] = 7.5
thresholds[capacities.index(13)] = 7.5
for i, t in enumerate(thresholds):
for j, g in enumerate(gate_improvements):
if g>t:
break
thresholdIdxs[i] = j
In [ ]:
for i in range(len(thresholds)):
thresholds[i] = 7.5
for i, t in enumerate(thresholds):
for j, g in enumerate(gate_improvements):
if g>t:
break
thresholdIdxs[i] = j
In [ ]:
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(distances)+1))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, axs1 = plt.subplots(int(len(capacities)/2),2)
axs1=axs1.flatten()
fig1.set_size_inches( 2*fig4.get_size_inches()[0],fig4.get_size_inches()[1]*len(capacities)/1.4)
for j, (c, ax1) in enumerate(zip(capacities, axs1)):
threshold = thresholds[j]
for i,d in enumerate(distances):
ax1.plot(gate_improvements, data["LogicalErrorRates"]["Forwarding"][d][j], linestyle='dashed', c=colors[i])
ax1.scatter(gate_improvements, data["LogicalErrorRates"]["Forwarding"][d][j], c=colors[i], label=f"d={d}")
ax1.vlines([threshold], ymin=[1e-7], ymax=[1e-0], colors=[colors[-1]], label="Threshold Point")
ax1.legend()
ax1.set_title(f'Threshold={threshold} for capacity={c}')
ax1.set_xlabel('Gate Improvement')
ax1.set_ylabel('Logical Error Rate')
ax1.set_xlim(1e-2, 1000)
ax1.set_ylim(1e-7, 1e-0)
ax1.set_xticks(gate_improvements)
ax1.loglog()
ax1.grid(which='major')
In [ ]:
def findMaxGi(_j, _d):
return min(
[
u
for u in range(len(gate_improvements))
if data["LogicalErrorRates"]["Forwarding"][_d][_j][u] == 0
],
default=len(gate_improvements),
)
cs = np.concatenate(
[
np.concatenate(
[
[c] * len(gate_improvements[thresholdIdxs[j] : findMaxGi(j, d)])
for d in distances
]
)
for j, c in enumerate(capacities)
]
)
ds = np.concatenate(
[
np.concatenate(
[
[d] * len(gate_improvements[thresholdIdxs[j] : findMaxGi(j, d)])
for d in distances
]
)
for j in range(len(capacities))
]
)
xs = np.concatenate(
[
np.concatenate(
[gate_improvements[thresholdIdxs[j] : findMaxGi(j, d)] for d in distances]
)
for j in range(len(capacities))
]
)
ys = np.concatenate(
[
np.concatenate(
[
[
data["LogicalErrorRates"]["Forwarding"][d][j][k + thresholdIdxs[j]]
for k in range(
len(gate_improvements[thresholdIdxs[j] : findMaxGi(j, d)])
)
]
for d in distances
]
)
for j in range(len(capacities))
]
)
xs = np.log(xs)
ys = np.log(ys)
import sklearn.linear_model
_model2 = sklearn.linear_model.LinearRegression()
def _model2paramsfunc(_cs, _ds, _xs):
_cs, _ds, _xs = np.array(_cs), np.array(_ds), np.array(_xs)
return np.column_stack([ np.where(_cs==c, _xs, 0) for c in capacities[:-1]]+[np.where(_ds==d,_xs,0) for d in distances]+[np.where(_cs==c, 1, 0) for c in capacities[:-1]]+[np.where(_ds==d,1,0) for d in distances[:-1]])
def _model2rootfunc(_cs, _ds, _desired_y):
_cs, _ds = np.array(_cs, dtype=int), np.array(_ds, dtype=int)
_numer_inter = np.zeros_like(_cs, dtype=float)+p0_hat
_denom_slope = np.zeros_like(_cs, dtype=float)
_i = 0
for _c in capacities[:-1]:
_denom_slope+=ps_hat[_i]*np.where(_cs==_c, 1,0)
_i+=1
for _d in distances:
_denom_slope+=ps_hat[_i]*np.where(_ds==_d,1,0)
_i+=1
for _c in capacities[:-1]:
_numer_inter+=ps_hat[_i]*np.where(_cs==_c, 1,0)
_i+=1
for _d in distances[:-1]:
_numer_inter+=ps_hat[_i]*np.where(_ds==_d,1,0)
_i+=1
return np.exp(
np.divide(
_desired_y - _numer_inter,
_denom_slope,
)
)
_model2.fit(_model2paramsfunc(cs, ds, xs), ys)
p0_hat, ( ps_hat) = (
_model2.intercept_,
_model2.coef_,
)
print( p0_hat, ps_hat)
ps_hat = ps_hat[2:]
_model2_errors_squared = (
ys
- _model2.predict(
_model2paramsfunc(cs, ds, xs)
)
) ** 2
print(np.sqrt(np.mean(_model2_errors_squared)))
x_bins = 100
y_bins = 100
_range = [[min(capacities), max(capacities)], [min(distances), max(distances)]]
counts, xedges, yedges = np.histogram2d(
cs, ds, bins=[x_bins, y_bins], range=_range
)
values, _, _ = np.histogram2d(
cs,
ds,
bins=[x_bins, y_bins],
weights=_model2_errors_squared,
range=_range,
)
counts = counts.T
values = values.T
with np.errstate(divide="ignore", invalid="ignore"):
values = values / counts
# scaled_to_error_values = (values/max(values.flatten()))*max(_model2_errors_squared)
scaled_to_error_values = values
scaled_to_error_values[counts == 0] = np.nan
x_centers = (xedges[:-1] + xedges[1:]) / 2
y_centers = (yedges[:-1] + yedges[1:]) / 2
X, Y = np.meshgrid(x_centers, y_centers)
points = np.column_stack([X.ravel(), Y.ravel()])
scaled_to_error_values_flat = scaled_to_error_values.ravel()
empty_bins = np.isnan(scaled_to_error_values_flat)
from scipy.interpolate import griddata
from scipy.spatial import cKDTree
if np.any(empty_bins):
interpolated_values_flat = griddata(
points[~empty_bins],
scaled_to_error_values_flat[~empty_bins],
points,
method="cubic",
)
empty_bins = np.isnan(interpolated_values_flat)
if np.any(np.isnan(interpolated_values_flat)):
tree = cKDTree(points[~empty_bins])
dist, idx = tree.query(points[np.isnan(interpolated_values_flat)])
interpolated_values_flat[np.isnan(interpolated_values_flat)] = (
interpolated_values_flat[~empty_bins][idx]
)
interpolated_values = interpolated_values_flat.reshape(scaled_to_error_values.shape)
def _model2scalefunc(_cs, _ds):
tree = cKDTree(points)
_, idx = tree.query(list(zip(_cs, _ds)))
return np.sqrt(interpolated_values_flat[idx])
print(np.mean(_model2scalefunc(cs, ds)))
In [ ]:
from scipy import stats
from scipy.optimize import minimize
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(distances)+1))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, axs1 = plt.subplots(int(len(capacities)/2),2)
axs1=axs1.flatten()
fig1.set_size_inches( 2*fig4.get_size_inches()[0],fig4.get_size_inches()[1]*len(capacities)/1.4)
fittings_dict = {c: {d: None for d in distances} for c in capacities}
for j, (c, ax1) in enumerate(zip(capacities, axs1)):
for i,d in enumerate(distances):
xs_for_cd = np.log(gate_improvements[thresholdIdxs[j]:findMaxGi(j,d)])
if len(xs_for_cd)==0:
continue
ys = np.log(data["LogicalErrorRates"]["Forwarding"][d][j][thresholdIdxs[j]:findMaxGi(j,d)])
res = stats.linregress(xs_for_cd, ys)
alpha_hat, beta_hat, stderr_s_hat, stderr_i_hat = res.slope, res.intercept, res.stderr, res.intercept_stderr
fittings_dict[c][d] = alpha_hat, beta_hat, (stderr_s_hat), (stderr_i_hat), np.sqrt(np.mean((ys-xs_for_cd*alpha_hat-beta_hat)**2))
def readout_stat(_rys_drawn, _xs_proj):
_res = stats.linregress(xs_for_cd, _rys_drawn)
_alpha_hat, _beta_hat = _res.slope, _res.intercept
return np.exp(_alpha_hat*(_xs_proj)+_beta_hat)
def ry_star():
return np.random.normal(loc=alpha_hat*xs_for_cd+beta_hat, scale=np.sqrt(stderr_s_hat*xs_for_cd+stderr_i_hat))
xs_projection = np.log(np.array([1e-2, 1000]))
readouts = np.array([readout_stat(ry_star(), xs_projection) for _ in range(10)]).T
los, his = np.quantile(readouts, [0.025,0.975], axis=1)
ax1.fill_between(np.exp(xs_projection), los, his, color=colors[i], interpolate=True,label=f"d={d}", linestyle="dashed", alpha=0.4)
ax1.plot(np.exp(xs_projection), np.exp(alpha_hat*xs_projection+beta_hat), color=colors[i])
ax1.scatter(gate_improvements, data["LogicalErrorRates"]["Forwarding"][d][j], c=colors[i], label=f"d={d}")
ax1.vlines([thresholds[j]], ymin=[1e-7], ymax=[1e-0], colors=[colors[-1]], label="Threshold Point")
ax1.legend()
ax1.set_title(f'Threshold={threshold} for capacity={c}')
ax1.set_xlabel('Gate Improvement')
ax1.set_ylabel('Logical Error Rate')
ax1.set_xlim(1e-2, 1000)
ax1.set_ylim(1e-7, 1e-0)
ax1.set_xticks(gate_improvements)
ax1.loglog()
ax1.grid(which='major')
In [ ]:
def rootfunc(_cs,_ds, _desired_y):
_rs = []
for c, d in zip(_cs, _ds):
if not fittings_dict[c][d]:
_rs.append(np.nan)
else:
_a_hat, _b_hat, _,_,_ = fittings_dict[c][d]
_rs.append(np.exp((_desired_y-_b_hat)/_a_hat))
return _rs
def _errorsSquared(_cs,_ds):
_rs = []
for c, d in zip(_cs, _ds):
if not fittings_dict[c][d]:
_rs.append(1.0)
else:
_, _, _, _, _mse = fittings_dict[c][d]
_rs.append(_mse)
return _rs
def _errorsFromFitting(_cs,_ds, _xs):
_rs = []
for c, d, x in zip(_cs, _ds,_xs):
if not fittings_dict[c][d]:
_rs.append(np.nan)
else:
_,_, _ss,_si,_ = fittings_dict[c][d]
_rs.append(np.sqrt(_ss*x+_si))
return _rs
def _predictFromFitting(_cs,_ds, _xs):
_rs = []
for c, d, x in zip(_cs, _ds,_xs):
if not fittings_dict[c][d]:
_rs.append(np.nan)
else:
_a_hat, _b_hat, _,_,_ = fittings_dict[c][d]
_rs.append(_a_hat*x+_b_hat)
return _rs
In [ ]:
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(distances)+1))
colors = [mplcolors.to_hex(c) for c in colors]
desired_y = np.log(1e-4)
fig1, ax1 = plt.subplots()
for i, d in enumerate(distances[1:]):
ax1.plot([k for k in range(len(capacities))], rootfunc(capacities, [d]*len(capacities), desired_y) , color=colors[i],label=f"d={d}")
ax1.scatter([k for k in range(len(capacities))], rootfunc(capacities, [d]*len(capacities), desired_y), color=colors[i])
ax1.legend()
ax1.set_title(f'Improvement Needed to Survive a 10,000 Rounds')
ax1.set_xlabel('Capacities')
ax1.set_ylabel('Gate Improvement')
ax1.set_yticks(gate_improvements)
ax1.set_xticks([k for k in range(len(capacities))])
ax1.set_xticklabels(capacities)
ax1.set_yscale('log')
ax1.grid(which='major')
In [ ]:
cs = np.concatenate(
[
np.concatenate(
[
[c] * len(gate_improvements[thresholdIdxs[j] : findMaxGi(j, d)])
for d in distances
]
)
for j, c in enumerate(capacities)
]
)
ds = np.concatenate(
[
np.concatenate(
[
[d] * len(gate_improvements[thresholdIdxs[j] : findMaxGi(j, d)])
for d in distances
]
)
for j in range(len(capacities))
]
)
In [ ]:
capacities_proj = np.linspace(2,60,100)
distances_proj = np.linspace(3,9,100)
cs_proj = np.array([capacities_proj for _ in range(len(distances_proj))])
ds_proj = np.array([[d]*len(capacities_proj) for d in distances_proj])
desired_y = np.log(1e-4)
gate_improv_proj = np.log(rootfunc(cs, ds, desired_y))
_range = [[min(capacities), max(capacities)], [min(distances), max(distances)]]
_counts, cedges, dedges = np.histogram2d(
cs, ds, bins=[100, 100], range=_range
)
_values, _, _ = np.histogram2d(
cs,
ds,
bins=[100, 100],
weights=gate_improv_proj,
range=_range,
)
_counts = _counts.T
_values = _values.T
with np.errstate(divide="ignore", invalid="ignore"):
scaled_to_improv_values = _values / _counts
scaled_to_improv_values[_counts == 0] = np.nan
c_centers = (cedges[:-1] + cedges[1:]) / 2
d_centers = (dedges[:-1] + dedges[1:]) / 2
C, D = np.meshgrid(c_centers, d_centers)
_points = np.column_stack([C.ravel(), D.ravel()])
scaled_to_improv_values_flat = scaled_to_improv_values.ravel()
from scipy.interpolate import griddata
_empty_bins = np.isnan(scaled_to_improv_values_flat)
if np.any(_empty_bins):
_interpolated_improv_values_flat = griddata(
_points[~_empty_bins],
scaled_to_improv_values_flat[~_empty_bins],
_points,
method="cubic",
)
_interpolated_improv_values = _interpolated_improv_values_flat.reshape(scaled_to_improv_values.shape)
fig, axs = plt.subplots(1,2)
fig.set_size_inches(20,10)
ax=axs[0]
contour_filled = ax.contourf(C, D, np.exp(_interpolated_improv_values), cmap='viridis', norm=mplcolors.LogNorm())
ax.contour(C, D,np.exp( _interpolated_improv_values), colors='black', norm=mplcolors.LogNorm())
cbar = fig.colorbar(contour_filled, ax=ax)
cbar.set_label('Gate Improvement')
ax.set_xlabel('Capacity')
ax.set_ylabel('Distance')
ax.set_yticks(distances)
ax.set_xlim(min(capacities),max(capacities[:-1]))
ax.set_ylim(min(distances),max(distances))
ax.set_title('Survive 10,000 Rounds')
x_bins = 100
y_bins = 100
_range = [[min(capacities), max(capacities)], [min(distances), max(distances)]]
counts, xedges, yedges = np.histogram2d(
cs, ds, bins=[x_bins, y_bins], range=_range
)
values, _, _ = np.histogram2d(
cs,
ds,
bins=[x_bins, y_bins],
weights=np.exp(_errorsSquared(cs, ds)),
range=_range,
)
counts = counts.T
values = values.T
with np.errstate(divide="ignore", invalid="ignore"):
values = values / counts
scaled_to_error_values = values
scaled_to_error_values[counts == 0] = np.nan
x_centers = (xedges[:-1] + xedges[1:]) / 2
y_centers = (yedges[:-1] + yedges[1:]) / 2
X, Y = np.meshgrid(x_centers, y_centers)
points = np.column_stack([X.ravel(), Y.ravel()])
scaled_to_error_values_flat = scaled_to_error_values.ravel()
empty_bins = np.isnan(scaled_to_error_values_flat)
from scipy.interpolate import griddata
from scipy.spatial import cKDTree
if np.any(empty_bins):
interpolated_values_flat = griddata(
points[~empty_bins],
scaled_to_error_values_flat[~empty_bins],
points,
method="cubic",
)
empty_bins = np.isnan(interpolated_values_flat)
if np.any(np.isnan(interpolated_values_flat)):
tree = cKDTree(points[~empty_bins])
dist, idx = tree.query(points[np.isnan(interpolated_values_flat)])
interpolated_values_flat[np.isnan(interpolated_values_flat)] = (
interpolated_values_flat[~empty_bins][idx]
)
interpolated_values = interpolated_values_flat.reshape(scaled_to_error_values.shape)
ax=axs[1]
contour_filled = ax.contourf( X, Y, np.exp(interpolated_values), cmap='viridis')
fig.colorbar(contour_filled, ax=ax, label='Squared Error')
ax.set_xlabel('Capacity')
ax.set_ylabel('Distance')
ax.set_yticks(distances)
ax.set_xlim(min(capacities),max(capacities[:-1]))
ax.set_ylim(min(distances),max(distances))
ax.set_title('Squared Error for Gate Improvement')
In [ ]:
_capacities = capacities[:-2]
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(distances)))
colors = [mplcolors.to_hex(c) for c in colors]
desired_y = np.log(1e-4)
fig1, ax1 = plt.subplots()
for j, d in enumerate(distances[1:]):
los, his = [], []
for i,c in enumerate(_capacities):
xs_for_cd = np.log(gate_improvements[thresholdIdxs[i]:findMaxGi(i,d)])
if len(xs_for_cd)==0:
los.append(0)
his.append(0)
continue
def readout_stat(rys_drawn):
_res = stats.linregress(xs_for_cd, rys_drawn)
_alpha_hat, _beta_hat = _res.slope, _res.intercept
return np.exp((desired_y-_beta_hat)/_alpha_hat)
def ry_star():
_cs, _ds = [c]*len(xs_for_cd), [d]*len(xs_for_cd)
return np.random.normal(loc=_predictFromFitting(_cs, _ds, xs_for_cd), scale=_errorsFromFitting(_cs, _ds, xs_for_cd))
readouts = np.array([readout_stat(ry_star()) for _ in range(10)])
lo, hi = np.quantile(readouts, [0.025,0.975])
los.append(lo)
his.append(hi)
ax1.plot([k for k in range(len(_capacities))], rootfunc(_capacities, [d]*len(_capacities), desired_y) , color=colors[j],label=f"d={d}")
ax1.scatter([k for k in range(len(_capacities))], rootfunc(_capacities, [d]*len(_capacities), desired_y), color=colors[j])
ax1.hlines(los+his, xmin=[k-1/2 for k in range(len(_capacities))]*2, xmax=[k+1/2 for k in range(len(_capacities))]*2, colors=[colors[j]]*(len(los)+len(his)), alpha=0.8)
ax1.vlines([k for k in range(len(_capacities))], ymin=los, ymax=his, colors=[colors[j]]*(len(los)+len(his)), alpha=0.8)
ax1.legend()
ax1.set_title(f'Improvement Needed to Survive 10,000 Rounds')
ax1.set_xlabel('Capacities')
ax1.set_ylabel('Gate Improvement')
ax1.set_xticks([k for k in range(len(_capacities))])
ax1.set_xticklabels(_capacities)
ax1.set_yscale('log')
ax1.grid(which='major')
In [ ]:
import scipy.stats
gi_offset = 0
_gate_improvements = gate_improvements[gi_offset:]
cmap = mpl.colormaps['plasma']
colors = cmap(np.linspace(0, 1, len(capacities)+1))
colors = [mplcolors.to_hex(c) for c in colors]
fig1, axs = plt.subplots(int(len(_gate_improvements)/2),2)
axs=axs.flatten()
fig1.set_size_inches( 2*fig4.get_size_inches()[0],fig4.get_size_inches()[1]*len(_gate_improvements)/1.4)
xs_projection = np.array([0, 50])
for j, (improv, ax) in enumerate(zip(_gate_improvements, axs)):
for i, c in enumerate(capacities):
if c==1 or c==30:
logicalErrors = np.array([data["LogicalErrorRates"]["Forwarding"][d][i][j+gi_offset] for d in distances])
ds = np.array(distances)
xs, ys = ds[logicalErrors>0], logicalErrors[logicalErrors>0]
ys = np.log(ys)
def logPr(params, ys_given, eps = 1e-10):
_alpha, _beta, _scale = params
_lik = stats.norm(loc=_alpha*(xs)+_beta, scale=_scale**2).pdf(ys_given)
return np.log(_lik+eps)
res = scipy.stats.linregress(xs, ys)
alpha_hat, beta_hat, scale_s_hat, scale_i_hat = res.slope, res.intercept, res.stderr, res.intercept_stderr
def readout_stat(rys_drawn, _xs):
res = scipy.stats.linregress(xs, rys_drawn)
_alpha_hat, _beta_hat = res.slope, res.intercept
return np.exp(_alpha_hat*(_xs)+_beta_hat)
def ry_star():
return np.random.normal(loc=alpha_hat*xs+beta_hat, scale=np.sqrt(((ys-alpha_hat*xs-beta_hat)**2)))
readouts = np.array([readout_stat(ry_star(), xs_projection) for _ in range(100)]).T
los, his = np.quantile(readouts, [0.025,0.975], axis=1)
ax.fill_between(xs_projection, los, his, color=colors[i], interpolate=True, label=f"c={c}", linestyle="dashed", alpha=0.4)
ax.plot(xs_projection,
np.exp(alpha_hat*xs_projection+beta_hat),
linestyle='--',
label=f"c={c}",
c=colors[i])
for j, (improv, ax) in enumerate(zip(_gate_improvements, axs)):
for i, c in enumerate(capacities):
if c==1 or c==30:
logicalErrors = np.array([data["LogicalErrorRates"]["Forwarding"][d][i][j+gi_offset] for d in distances])
ds = np.array(distances)
xs, ys = ds[logicalErrors>0], logicalErrors[logicalErrors>0]
ax.scatter(xs, ys, c=colors[i])
ax.legend()
ax.set_title(f'Projection for {improv}x Gate Improvement')
ax.set_xlabel('Distance')
ax.set_ylabel('Logical Error Rate')
ax.set_ylim(1e-8, 1e2)
ax.set_xlim(*xs_projection)
ax.set_yscale('log')
ax.grid(which='major')